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ABSTRACT 

The correct qualitative features of the helium n elementary excitation spectrum are 
derived microscopically for a realistic interatomic potential. The strong repulsive core 
is included by using a reaction matrix in the Hamiltonian. The attractive well is success- 
fully included by assuming a generalized Bose-Einstein condensation. The pair Hamilton- 
ian is diagonalized by the thermodynamically equivalent Hamiltonian method. Numerical 
solutions yield spectra with phonon and roton regions. The spectrum energies are too 
high for the Yntema- Schneider potential. Another potential, constructed to fit virial co- 
efficient data classically, gives better results. Spectra are presented for a series of at- 
tractive well strengths. 
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THE EXCITATION SPECTRUM FOR A BOSE GAS WITH 
REPULSIVE AND ATTRACTIVE INTERACTIONS 


by Gerald V. Brown 
Lewis Research Center 

SUMMARY 

The correct qualitative features of the helium II elementary excitation spectrum are 
derived microscopically for a realistic interatomic potential. The strong repulsive core 
is included by using a reaction matrix in the Hamiltonian. The attractive well is success- 
fully included by assuming a generalized Bose-Einstein condensation. The pair Hamilton- 
ian is diagonalized by the thermodynamically equivalent Hamiltonian method. Numerical 
solutions yield spectra with phonon and roton regions. The spectrum energies are too 
high for the Yntema- Schneider potential. Another potential, constructed to fit virial co- 
efficient data classically, gives better results. Spectra are presented for a series of at- 
tractive well strengths. 


INTRODUCTION 

The attempts which have been made to explain the unusual properties of superfluid 
helium (He II) have met with success or with difficulties, depending on the level of the ap- 
proach. On the successful side are the phenomenological theories with empirically ad- 
justed parameters. On the less successful side are the microscopic theories which at- 
tempt to derive the properties of liquid helium II from the interatomic forces of helium 
atoms. No microscopic theory has been able to deduce the energy spectrum (energy as a 
function of momentum) for a realistic potential. Prior to the present work, not even the 
qualitative features of the spectrum (the phonon-like part and the roton minimum) have 
been successfully derived from an interatomic potential including both strong repulsion 
and an attractive well. 

The phenomenological derivations of the spectrum are well known. Landau (refs. 1 
and 2) deduced from specific heat data that the spectrum of elementary excitations in he- 
lium should contain a phonon- like part (a linear portion at and near zero momentum) and 



a group of higher energy excitations, which he termed ''rotons. " The energies of these 
excitations were postulated (refs. 1 and 2) to be 
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where p is momentum and c, p', ji, and A are constants adjusted to fit the experimen- 
tal specific heat data. The phonon- like excitations contribute a T term to the specific 
heat, whereas the roton excitations make an exponential contribution because of the "en- 
ergy gap" and because Boltzmann statistics is satisfactory for the rotons. The previous 
relations (eq. (1)), with the constants appropriately adjusted, give very good qualitative 
agreement with the data from neutron scattering experiments (refs. 3 to 5) performed 

years after Landau's papers. Feynman (ref. 6) derived a similar energy spectrum from 

4 

more basic principles, with an argument based on the Bose-Einstein statistics of He 

atoms. Because Landau's derivation took no account of statistics, it made no qualitative 

4 3 

distinction between He and He' . Feynman's work is not entirely microscopic, however, 
for he utilizes the experimentally determined structure factor (ref. 7) for the liquid. 

Two groups of experimenters (refs. 3 to 5) have measured the spectrum of elemen- 
tary excitations in helium by neutron scattering experiments as proposed by Cohen and 
Feynman (ref. 8). The neutrons are scattered by density waves in the liquid. According 
to Pines (ref. 9), these density waves have the same energy spectrum as the elementary 
excitations for a system of bosons. The spectrum is shown in figure 1, and the phonon 



Figure l. - The experimental spectrum. 
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and roton regions are identified there. 

A chief aim of the many-body problem for helium is to derive this spectrum from the 
theoretically or experimentally determined interatomic forces of helium atoms, without 
further input from experiments. This work derives the spectrum from a realistic inter- 
particle potential by combining three methods. First, a generalized or ’’smeared” Bose- 
Einstein condensation is assumed because the scattering length of the interparticle poten- 
tial is negative. Second, a modified reaction matrix is introduced to handle the strong 
repulsive core of the helium potential. Third, a Thermodynamically Equivalent Hamilto- 
nian (TEH) method permits the inclusion of all forward, exchange, and pair scattering in- 
teraction terms. 


THE IMPERFECT BOSON GAS 

The microscropic problem of the imperfect boson gas has been considered by numer- 
ous authors (refs. 10 to 23). The starting point is the second-quantized Hamiltonian for 
a system of bosons having an interparticle potential operator V 

Z_-/ k l k 2 k 3 k 4 

(Symbols are defined in appendix A. ) All indices are understood to be vectors, although 
the vector signs have been suppressed to avoid overcrowding. The system is enclosed in 
a box of volume ft, and the summations run over all allowed free-particle states in the 
box. The operator aj£ is the creation operator for a plane wave state with propagation 
vector k, and a^ is the corresponding destruction operator. The interaction potential 
v(r) is a spherically symmetric function of the distance r between two atoms, and ap- 
pears here in matrix elements with respect to free-particle two-body states, for example, 
(kjkg | v|kgk 4 ) . These elements can be expressed in terms of matrix elements with re- 
spect to one-particle states in a central potential by changing to center-of-mass coordi- 
nates. The result is 

<k,k 2 |v|k 3 k 4 ) = 6 kj+k 2 , k 3 +k/^3 - k l> = 6 k 1+ k 2 , k 3+ k/<1> 

where q = |kg - kjJ is the momentum transfer, 6 is the Kronecker delta, and v(q) is 
1/ft times the three-dimensional Fourier transform of the real- space potential 
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v(q) = — f d^r e^’ r v(r) 

n J 

In the present treatment, the change to center-of-mass coordinates must be postponed 
until after the K matrix is introduced. The form (kjkgfvlkgk^) will be retained until 
that point. 

The diagonalization of the Hamiltonian (eq. (2)) has not been accomplished. Many 
authors drop most of the interaction terms, although attempting to keep as many as pos- 
sible. These retained terms should be diagonalizable or amenable to some other treat- 
ment. One successful theory based on a "truncated" Hamiltonian is the Bardeen, 

Cooper, and Schrieffer theory of superconductivity (refs. 24 and 25). 

The Hamiltonian to be used here is the "pair Hamiltonian" of Girardeau and Arnowitt 
(ref. 21). All interaction terms are neglected except the forward scattering, exchange 
scattering, and pair scattering terms. These three types of terms are represented in 
figure 2 and have the following second-quantized forms, (qp| V|qp) a^a-pa^ap, 
(qp|v|pq}a+apapaq, and (-qq|v|-pp)a^a*q_pap, respectively. The truncated Hamilto- 
nian, now called the pair Hamiltonian Hp (ref. 21), has become 


H 


P 



£ <®l V lw> a q a pV p 
qp 


+ i Yj < qp l V l p< D a q a ^p a q + ;£ ( “QQ I V | - pp) a ^q a q ; 
qp qp 


~ a a 
■a a -d d 
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Particles with wave numbers as labeled 

Interaction labeled by appropriate matrix element 



Forward scattering Exchange scattering 
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particles 

\<-qq|V|pp)l 


-PI , .IP 

Incoming 
particles 

Pair scattering 


Figure 2. - Interaction terms retained in pair Hamiltonian. 
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The restrictions on the sums are necessary to prevent duplication of terms. 

Even the simplified Hamiltonian (eq. (3)) has not been diagonalized. If two additional 
simplications are made, diagonalization can be achieved: These are: (1) keep only those 
interaction terms containing at least two creation or annihilation operators with subscript 
zero; (2) approximate both a* and a Q by ^N. If the potential is repulsive, that is, 
v(0) > 0, the resulting "Bogoliubov Hamiltonian” can be diagonalized (ref. 10). The sec- 
ond simplication, called the Bogoliubov approximation, is justified for weak interactions 
near absolute zero because nearly all particles are expected to be in the zero-momentum 
state. There is a canonical transformation on the single-particle operators which diago- 
nalizes the Bogoliubov Hamiltonian. Such a transformation is called the Bogoliubov 
transformation and will be used later in the present work. A phonon spectrum for low 
momenta occurs in this approximation and a free-particle spectrum at high momenta. In 
the intermediate range of momenta the spectrum has a region connecting the linear and 
quadratic sections. For an appropriate repulsive potential this region could have the gen- 
eral shape of the roton region of the liquid helium spectrum. The model thus has two fea- 
tures resembling liquid helium: the low-momentum phonon spectrum and at least a hint 
of a roton region. 

Girardeau and Arnowitt (ref. 10), using a variational method, consider the entire 
pair Hamiltonian (eq. (3)), without making the Bogoliubov approximation. Wentzel 
(ref. 14) and Luban (ref. 15), also studying the Hamiltonian (eq. (3)), allow thermal ex- 
citation, and find a simpler diagonalizable Hamiltonian which gives the same thermody- 
namics as that of equation (3). All three of these studies find an energy gap in the low- 

momentum excitation spectrum. That is, E(0) = 0, but lim E(k) ¥= 0. These results are 

k— 0 

for weak interparticle potentials with v(0) > 0, and they assume that Bose-Einstein con- 
densation takes place with particles "condensing” into the zero-momentum state. The 
spectrum of this pair Hamiltonian model, which includes more terms than the Bogoliubov 
Hamiltonian, is, nevertheless, further from that of helium II. 

There are two important differences between the interparticle potential used in the 
aforementioned studies and the actual helium potential. First, the helium potential is 
much more strongly repulsive at close approach. Secondly, the scattering length for he- 
lium is negative. (In fact, the scattering length is so negative that the atoms can almost 
form a two-body bound state. ) For the Yntema-Schneider (Y. S. ) potential (refs. 26 and 
27), the scattering length is negative, and its magnitude is several times the repulsive 
core diameter. For a potential with a singular repulsive core, as in helium, the sign of 
the scattering length plays the role that the sign of v(0) plays in a weak potential, roughly 

*The scattering length characterizes the "net effect” of a potential, in that a positive 
scattering length indicates net repulsion, and negative scattering length indicates net at- 
traction. The scattering length for hard spheres is equal to the diameter of a sphere. 
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speaking. Thus a system of atoms with negative scattering length may be expected to 
correspond more nearly to a system with v(0) < 0 than to one with v(0) > 0. The repul- 
sive core and the negative scattering length each require changes from the methods used 
by the authors mentioned in the preceding paragraph. 


THEORETICAL TECHNIQUES 

Three main methods will now be described, which in combination allow a potential 
with singular repulsive core and an attractive well to be studied using the pair Hamilto- 
nian. To deal with potentials with negative scattering length (like helium), a generalized 
or "smeared” Bose-Einstein condensation is assumed. To avoid the infinite matrix ele- 
ments of the repulsive core, a type of reaction matrix is used in place of the potential. 
Finally to obtain a diagonalized form from the pair Hamiltonian, the Wentzel thermody- 
namically equivalent Hamiltonian (TEH) method is used. The three methods are dis- 
cussed in the three following sections. 


Generalized Bose-Einstein Condensation 

Consider first the effect of having a potential with v(0) < 0. Girardeau (refs. 22 and 
23) argues that for weak attractive potentials, where v(0) < 0, the lowest energy state is 
not one in which the zero-momentum state contains a finite fraction of the particles 
(simple Bose-Einstein condensation). It is rather one in which a large number of distinct 
low-momentum states contain a finite fraction of the particles but any single state con- 
tains only a negligibly small fraction. In this generalized condensation all states with 
momentum less than a cutoff momentum p Q are assumed to have zero energy. A finite 
fraction of the particles occupy the group of states, but no single state, not even the 
p = 0 state, contains a finite fraction. The group of states containing condensed parti- 
cles draws arbitrarily close to zero momentum and bears a deceptive resemblance to 
simple condensation. It is the fact that no single state is macroscopically occupied (that 
is, contains a finite fraction of the total particles) that distinguishes the two types of con- 
densation. 

In reference 17, the descriptive statements about generalized condensation are for- 
mulated as follows. As the thermodynamic limit is taken, that is, N — °° but N/fi is 
constant. 


— oo 
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P o = 0(fi(- z >) (z > 0) 

( a p a p> = 0(n z ^) for I p I < P 0 and 0 < z (p) = 0(1) < 1 

and 

I <a£y = o<n> 

IXPq 

The results of three studies support the use of a smeared condensation. Girardeau 
demonstrates by variational means that smearing gives a lower energy for weakly attrac- 
tive systems than does a simple condensation. Sawada and Vasudevan (ref. 28) show in a 
simplified model with negative scattering length that the states into which the particles 
condense should be a combination of zero- and nonzero-momentum states. This combi- 
nation of states may be viewed as a smeared condensation. 

Additional incentive to investigate the smeared type of condensation is provided by 
the work of Luban (appendix E of ref. 17). He showed that in the pair Hamiltonian model 
with a hard core pseudopotential and weak attractive interactions, a smeared condensa- 
tion leads to a phonon-like spectrum for low-momentum excitations rather than to the 
energy gap predicted by simple condensation. The modifications to be made below to 
treat strong potentials do not change the character of these low-momentum excitations if 
the scattering length is negative. Thus it appears reasonable to use a smeared condensa- 
tion in a study of helium, which has a negative scattering length. 


The Reaction Matrix 

To treat an interparticle potential with a strong repulsive core, the pair Hamiltonian 
(eq. (3)) is inadequate as it stands. The matrix elements of the interaction potential that 
appear in equation (3) are arbitrarily large for arbitrarily strong cores. This is easily 
seen by considering a repulsive core of uniform height V Q and radius a as shown in 
figure 3. Then 

4ttV /*a 

(kj + q, k 2 - qlvlk^) = v(q) =— — r sin (qr) dr 
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a 

Separation, r 


Figure 3. - Square repulsive 
core. 


This general matrix element is proportional to V . In the limit of V Q — 00 (hard core), 
progress can still be made by summing enough terms in a many-body perturbation series 
of the exact Hamiltonian to obtain a finite result. 

The method used here to effect this summing is similar to that of Brueckner and 
Sawada (refs. 14 and 15). The matrix elements of the interparticle potential operator are 
replaced by the elements of a type of reaction matrix. Since the reaction matrix K is 
defined by 


K = V + VGK (4) 

where V is the exact two-body potential operator and G is a Green's function operator, 
it has an expansion of the form 

K = V + VGV + VGVGV + . . . 

The use of matrix elements of K in place of matrix elements of V (which is the first 
term in the expansion of K) brings many more interaction terms into the Hamiltonian 
without complicating its form. If the matrix elements of K are calculable from V, then 
the use of K elements effectively presums enough interaction terms to give a finite re- 
sult. (Note that, for very weak V, the K operator approaches V but that, for singular 
core potentials, elements of K with respect to plane waves are still finite whereas those 
of V are infinite (refs. 12 and 13).) 

The arguments for this replacement of V by K and the selection of the operator G 
are made in appendix B by considering the many-body perturbation expansion of the free 
energy. The perturbation expansion of the free energy based on the pair Hamiltonian con- 
taining K's instead of V's is more nearly like the expansion of the free energy based on 
a complete nontruncated Hamiltonian. The operator G must be appropriately chosen, 
however. The argument in appendix B shows that G = -l/H Q , where H Q is the kinetic 
energy operator, can be used. 
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Matrix elements of K with respect to two-body plane-wave states are needed to in- 
sert into the truncated Hamiltonian. To reduce the calculation of the matrix elements to 
manageable proportions, an approximation is made - the center-of-mass approximation. 
The matrix elements of K with respect to two-body states are approximated by elements 
with respect to one-body center-of-mass states. Details of the center-of-mass approxi- 
mation, the types of matrix elements needed, and the decomposition of the elements into 
partial waves are contained in appendixes C and D. The important results from appen- 
dixes B, C, and D are the following: 

(1) Matrix elements (kp|v|qr) in the Hamiltonian are to be replaced by the corre- 
sponding reaction matrix elements (kp|K|qr) . 

(2) The two-body elements (kp|K|qr) are to be approximated by the one-body center- 
of-mass elements 



The integral equation K = V + VGK for these one-body elements is decomposed into a 
similar equation for each partial wave by expanding all elements of the matrices in 
spherical harmonics. The integral equation for each partial wave can be solved by ma- 
chine. Only even partial waves are needed, and three of these give sufficient accuracy 
in the energy spectrum from zero momentum to just past the roton minimum. 


The Thermodynamically Equivalent Hamiltonian 

The modification of the truncated Hamiltonian (eq. (3)) to allow treatment of strongly 
repulsive cores has not changed its basic form. It is still of the type which can be 
treated by the TEH method (refs. 13 and 14). In the first part of this section the meaning 
of "thermodynamically equivalent" is discussed, and especially the question of how this 
method can be applied in the present work where the use of the reaction matrix is justi- 
fied only in the limit as T — 0 (appendix B) . 

The essence of the Wentzel method is that a simpler Hamiltonian than equation (3) 
can be found which gives the same partition function as equation (3) in the limit as 

S2 — °° (the thermodynamic limit) but which can be exactly diagonalized by the Bogoliubov 
o 

transformation. Two systems with the same partition function have exactly the same 
thermodynamics, but, in general, this does not guarantee any microscopic similarity. 
The microscopic similarity is a central point of this work, however, which attempt to 
o 

Note that the Bogoliubov transformation and the Bogoliubov approximation are dis- 
tinct. 
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find the energy spectrum of elementary excitations (normal modes) in a dense boson gas 
resembling helium. The relation, if any, between the spectrum of excitations found by 
diagonalizing the thermodynamically equivalent Hamiltonian and the experimentally deter- 
mined spectrum of excitations (refs. 3 to 5) for liquid helium must be examined. 

The grand partition function of a general interacting system is 

'^state-^o 
2(J3, p) = \ e state ° 

all dis- 
tinguish- 
able states 

where is the total number operator. If the system Hamiltonian (eq. (3)) could be 

diagonalized so it could be written as 


H = £ E k“J“k 

k 


where a* creates |k), a state with momentum k and energy E^, then the grand parti- 
tion function would be 


Wm) 



( 5 ) 


where ^ indicates a sum over all possible sets of occupation numbers n^. The en- 

{n k } 

ergies in this last expression are temperature independent numbers. The Wentz el 
method diagonalizes the truncated Hamiltonian (eq. (3)) in the sense that 
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( 6 ) 


Here, as in equation (5), the summation runs over all sets of occupation numbers n^. 
This partition function (eq. (6)) is constructed equal to that of equation (5) for any ( 3 , but 
in general e^(/3), a function of temperature, bears no simple relation to of equa- 
tion (5). If, however, for a range of temperature from zero to some finite temperature 
Tf, is independent of temperature, then it is expected that 

£ k<« = E k - ^ ' € k 

for T < Tf. 

The integral equations that determine 6^(0) are not very sensitive to temperature 
near T = 0 because temperature enters ohly in the thermal expectation values of the 
number operator and of the 

limiting and, in general, nonzero values as T — 0. For very small T, departures from 
the limiting values will be arbitrarily small. The reason for this can be most easily 
seen after the solutions of the equations are obtained. The lowest energy excitations are 
seen to be phonon-like, that is, their energies are proportional to momentum. This type 
of excitation spectrum is much more "rigid” against thermal excitation than, for ex- 
ample, a quadratic free-particle spectrum. In the former case many fewer states have 
energies of the order of kT than in the latter case. The relative rigidity against ex- 
citations causes the number of excited particles to be relatively constant near zero tem- 
perature. This, in turn, leads to the insensitivity to temperature of the excitation spec- 
trum and all thermal expectation values as T — 0. Thus the phonon- like spectrum leads 
to the temperature independence of e^O?) as /3 — 00 (T — 0) and hence to the assertion 
that the e^OS — °°) found by Wentzel’s method is the same as the normal mode spec- 
trum. 

As previously mentioned, the experimentally determined excitation spectrum is not 
strictly temperature independent in the temperature range 1. 1 to 1. 8 K. The roton min- 
imum is 5 percent lower in energy at 1. 8 K than at 1. 1 K. It is therefore questionable 
whether the experiment measures purely normal (noninteracting) modes, or perhaps 
whether strictly normal modes even exist. However, the temperature dependence is not 
strong, and at least part of it may be due to the slight change in density, so the modes 
measured experimentally are, at worst, weakly interacting. 


pair destruction operator ( a ^ a _j J ) • These will have 
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INTEGRAL EQUATIONS FOR THE SPECTRUM 

Consider again the pair Hamiltonian (eq. (3)), replacing the V-matrix elements with 
the corresponding K-matrix elements to yield the following: 


X, ' ") “k"* + 2 Z <<JPl K I<lP> a J a >q a j 


+ I £ (pq|K|qp)a^a ii a p+ I^ (-„|K|-pp)a2 q a+a_ p a f 


fi 2 k 2 


\ ) a k a k + 2 E K qpqp a q a q a p a p 


TT" "F *F 

K ala aTa„ 


+ i V k a + a + a a + - V K 
% A pqqp p q q p o Li - 


-f- 4 * 

a a a a 
qq-PP -q q -p P 


where = (pq|K|rs> . (Note that ± £ K qpqp a q a p a q a p = \ E K qpqp a qVp a p " 

pq pq 

1 V K kkkk a^a k . \ The form (eq. (7)) can be treated by the TEH method (ref. 13). As ii 

2 k 

Luban’s version (ref. 14) of Wentzel’s method, \i, the chemical potential, has been in- 
serted in Hp. This is soon eliminated from the equations. 

To find the TEH, first define new operators BjJ, C k , and C+ in the following 


~ ®k + ^k 


a -k a k “ C k + ^k 


ii ii i n i i ■ Mil 1 1 1 ii ail ii ■ in 
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where the values of the real c-numbers ^ and 77^ will be chosen to eliminate some 
terms from the new equivalent Hamiltonian. These substitutions give 




2m 2 / K K 2 


Vhs + o X/ K qpqp^ B p + V^ B q + V 


pq 


+ , Z K pqqp^ B q + ^q^ B P + V + O 2 K -qq-PP^ C q + V^ C P + V P* 


pSq 


pq 

p*q 


Because = B^, 


h p = 


— — - 11 - 
2 m 


K 


kkkk\ +„ . 1 


aja, +- V K iB^B^ + f (a^a - £ ) 
kTc 9 Lu qpqi> P q Tqq s q' 


pq 


+ ^ a p a p " V^q + ^p^q^ % E K pqqp^- B q B p + ^ a p a p " V + ^ a q a q ~ V^p + ^p^ 


pq 

p*±q 




K 


qq-ppI C q C 


_ + 77 „(a a _ 

p 'q v p -p 


V 


+ (a+a+ 

v q -q 


V T? P + Vq] 


pq 

p^±q 


Note K = K and K = K ; that is, the reaction matrix element is sym- 
qpqp .oq qppq pqqp' 

metric with respiV to interchange of the first and second pairs of indices. Using these 
facts and manipulating the dummy indices yield 
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fi 2 k 2 

2m 


K, 


M 


kkkk l 


is 

pq 


K ,(BpBg + 2 {„a„\ - {„£„) 


qpqp 


’pqq ’p’q' 


2 2 K qPPq^ B q B P 2 ^P q q 
pq 

P^±q 

+ i Z K -qq-PP [C q C P + V^-q + " Vp 1 

^ pq 
p^q 

Now collect the terms into two groups, putting those containing B, B + , C, and C + 

into 


H i 9 Z K qpqp B p B q + 9 Z K qppq B q B p + 9 Z K -qq-pp c q c i 


qp 


pq 

p^±q 


pq 

P*±q 


All the other terms are put into 


H TE ~ u + 



K 2 k 2 _ K kkkk\„+ 
2m 2 


/ qp 


+ Z K qppq^p a q a q + o Z K -qq-pp r? p^ a q a -q + a q a -q' ) 


qp 

q*±p 


qp 

q^p 


where 


u " 9 Z K qpqpVq " 9 Z K qppq^p^q " 9 Z K -qq-PiWq 


qp 


pq 

P^±q 


qp 

q*p 
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which is a c-number. Let 


f k = 


R 2 k 2 

2m 


K, 


M “ 


kkkk 


+ ^kpkp^p + 


Z 

p#±k 


K 


kppk^pl 


and 


(8) 


h k ^ K -kk-pp 7? p 
p*k 


Then 


«te - v * Z + \ Z h k ( < a -k + V-k> (9) 

k 1 k 

According to Wentz el’s TEH method, H TE will lead to the same thermodynamic 
properties as Hp if i- and are chosen as £p = ( a p a p) an< * ^p = ^ a p a -p^ = ^ a -p a p^ 
where the bracket denotes an average with respect to the grand ensemble. It has been 
shown that the Wentz el result remains valid for the "Hamiltonian” which has resulted 
from replacement of V-matrix elements by K-matrix elements. 

In appendix E the Bogoliubov transformation is used to diagonalize equation (9), the 
result being 


h TE - U o - Z 6 k“k“k 

k 


where 


and 



e k = 



( 10 ) 


(ID 
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It is noted in appendix E that the Bogoliubov transformation is not valid for any k for 
which e k = 0. In a smeared condensation there is a group of states with = 0, but 
these states are grouped arbitrarily close to k = 0 in the limit of infinite system volume. 
That is, these states all have wave numbers less than an arbitrarily small cutoff wave 
number p Q . All the equations are to be solved only in the limit of infinite volume. Thus 
the Bogoliubov transformation can be performed for all states not in the condensate, that 
is, that have finite momenta in the volume limit. The generalized or smeared Bose- 
Einstein condensate exists in states with e = 0, and it is not necessary to perform any 
transformation to find the energy of these states. The number of particles in these states 
is found by taking the difference between the number residing in excited states and the 
total number. 


In equations (8) the limit as k approaches zero gives (denoting lim f(k) as f(0) and 

v k— 0 


lim h(k) as h(0)) 
k— 0 7 


f <°> = - i K oooo + I K opop^p + Z K oppo« 

* P p*0 


h <°> = E ^o-pp’jp 

p#0 


(12) 


Because e(0) = Vf 2 (0) - h 2 (0) = 0, f(0) = ±h(0). Paralleling Luban (ref. 15), let f(0) = 
-h(0). Then 71 can be eliminated from equations (8) using equations (12) 


9 K oooo “ X K opop^p " X K oppo^p " X ^oo- pp^ i 


P 

p^O 


P 

p^O 


— K - K 
2 oooo oooo 


E v E <w 


p<p„ p<p c 

p^O 


X K opop^] 


P>P r 


V K ^ - V K 77 
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Then equations (8) become 


f k “ 


K 2 k 2 

2m 


+ Z V K koko 

P<Po 


+ K kook “ K oooo) ~ ^oooo X/ + 

P<P 0 

p^O 


+ X/ ^ K pkpk + K pkkp " K popo ” K poop^p " X/ ^oop-p^p 
P>P C P>P 0 


^k ^ ^k-koo^p + ^k-kp-p^p 


P<P r 


P>P C 

p*k 


where terms with p < p Q are separated from the sums, and continuity of the K-matrix 
elements near zero momenta with respect to any of the indices is used. Two terms, 
K 00 oo / 2 anc * ”^kkkk/ 2 ’ ^ ave t* een dropped from f^ because they cancel in the center 
of- mass approximation as can be seen from appendix C. 

It was shown by Luban (ref. 15) that (in the thermodynamic limit) as p — 0, 

T? p ~Sp + -. Then since | p = (aja p ) , the sums £ p and ^ 


T/ p are merely the 


P<P 0 P<P 0 

grand ensemble averages of the number of particles in the "smeared” condensate. If 
this number is N Q and the sums are replaced by integrals by letting 



then 


+ NoOWo + »W - 3K oooo> + — 


2m 


(2 v)' 


/ ( [ k p 


ikpk + ^pkkp " ^popo ” ^poop^p 



h, = N K, , + • 

k o k-koo 


n 




J" ^k-kp-p^p ^ p 


Note that in the integral form it is permissible to ignore the restrictions p * k, p + ±k, 
and p > p 0 and to carry out the integrals throughout p-space. The first two restrictions 
do not affect the integration since the excluded discrete states, p = k or p = ik, make 
negligible contributions. In the thermodynamic limit p Q is assumed to approach zero, 
and the interval of integration that would be excluded by the last restriction makes a neg- 
ligible contribution. 

From appendix E, 



— coth 
e k 



-1 


and T? k 



The one remaining relation needed is N = ^ ( a^a^ ) = ^ £ k or 

k k 

N q = N - [fi/(27r) 3 ] J d 3 k. 3 

These equations constitute a set of nonlinear coupled integral equations. Once the 
necessary elements of the reaction matrix K are calculated, the equations may be 
solved numerically for e(k), the quasiparticle energy spectrum, and for N Q /N, the 
fraction of particles (not quasiparticles) in the condensate. 

It is convenient for machine solution (and for simplicity of form) to nondimensional- 
ize the integral equations. This can be done by making the following definitions: 


x^ka y = pa Y = ~ 

a 


F(x) . l™? f H (x) a 2ma! h E(x) = Ima 2 


ft 


ft 


ft 


P a - a 3 P sUsa 3 B a — 

n n 2ma 2kT 


3 Here it is implicitly intended that, for p < p Q , £ p and ?] p inside integral signs 
mean simply the smooth extrapolations from p > p Q . The actual values of £ p and 77 p 
for p < p Q are much larger and their contributions have been separated out already. 
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ma 2 



v(r) 


° —in 

where a is any reference dimension, taken to be 2. 3 A (2. 3x10 m), the approximate 

core size (fig. 16), in this work. In these dimensionless quantities, the equations are 


F(x) x + 87rP o (jf xoxo + Jf XQOX - 3j ^ 000 ) 


hf' 


+ 2 J d y ^y^ J *xyxy + J ^cyyx _ **yoyo _ ^ooy^ 

7 T 


- 1 f 

IT • / 0 


y dyJf ooy-y\ 


H(x) = 87rP o Jf x _ xoo + i ^ y 2 dy * x _ xy _ y ^ 


r ( ls ) 


E(x) = VF 2 (x) - H 2 (x) 


p g = p - -- r * 

° 2 7t 2 y ° 


,y 2 dy 


As discussed in appendix C a general element K X y ZU of the reaction matrix is approxi- 
mated in this work by the one-particle, center-of-mass reaction matrix element 
K x _ y z _ u . Hence, in calculations K xoxo is replaced by x , for example. The equa- 

2 2 2 2 
tions for the one-particle reaction matrix elements, derived in appendix D, are these 


<k|K|p> = J] {21 + 1)P^(£ • p) <k|K|p) z 
l 


where 


o z * 00 

<k|jr|p = <k|f"|p> z -ZJ d(qa) (k|^|q) z (q|jr|p> z 


(D6) 


and 
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yV(r)i z (kay)j z (pay) dy 


Equations (13) are coupled and nonlinear but can be solved by relatively simple nu- 
merical methods on a computer. The same basic method is used as for solving the reac- 
tion matrix integral equation. Although the equations for F, H, and E are coupled and 
nonlinear, their solution is less demanding of machine computing time and memory stor- 
age than the solution of the linear reaction matrix equation (eq. (D6)). One reason is 
that F, H, and E are only one- dimensional arrays as compared to the two-dimensional 
jt and require two orders of magnitude less storage. Secondly, the integrands in the 
F, H, and E equations vanish much more rapidly with large momentum and the numer- 
ical integration may be stopped sooner. 

The method of solution is to start with a reasonable guess for the functions F(ka), 
H(ka), E(ka), and P Q = P and to insert these quantities in the right-hand members of 
equations (13). The new values for the functions calculated by doing the integrals were 
then averaged with the original guesses, the result being used as the next approximation. 
With a reasonable initial guess, this iterative method converged to within a fraction of 
1 percent in about ten iterations. The reason for averaging new values with old ones was 
to prevent oscillations around the actual solution. 

Because all the integrands decrease very rapidly after the roton minimum is passed, 
it is necessary to carry the integration only to p/ft ~ 4 A 4 (4x10^ m *). (The roton 
minimum is observed experimentally at p/ft = 1.8 A 4 or 1. 8x10^ m . ) 


THE INTERPARTICLE POTENTIAL 

To solve the integral equations numerically, a specific potential function v(r) must 
be chosen. The results presented in the next section will show great sensitivity to the 
strength of the potential well. Unfortunately, the well region has not been determined 
accurately by either experiment or theory. Consequently two potentials are presented in 
this section. They both fit the measured virial coefficients. One is an accepted poten- 
tial; the other is constructed simply to illustrate the effect of a weaker attractive well. 

Probably the best known expressions for the helium potential are the Slater- 
Kirkwood potential (ref. 29) 

v(r) = (770e -4 ‘ 6r - 1. 49 r" 6 )xl(T 12 (14) 

and the Yntema- Schneider (Y. S. ) potential (ref. 27) 
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v(r) = (1200e 


(15) 


-4. 72r _ j 24 r -6 _ j 89 r -8 )xl0 12 

where v(r) is in ergs and r is in angstroms. The Y. S. potential is shown in figure 4. 

The former was derived on theoretical grounds. The attractive tail is calculated 
from second-order perturbation theory for the Van der Waals interaction of two neutral 
atoms. The form of the repulsive core, due to Slater (ref. 30), is a first approximation 
to the overlap energy of two atoms which are close together. The attractive part domi- 
nates for large interparticle separation and the repulsive part for very small separation. 
The potential for intermediate separations, in particular in the vicinity of the minimum 
of the potential well, is not determined with any great accuracy, but in fact is the result 
of adding the repulsive and attractive terms together in the range of intermediate separa- 
tion. 
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The Y. S. potential was derived from experimental measurements of the second virial 
coefficient between 273 and 1473 K (ref. 26). The form ae~^ r - cr"® - dr”^ was as- 
sumed. The value for c was taken from a theoretical derivation of London (ref. 31) and 
the value of d from Margenau (ref. 32), who calculated it to correspond to the London 

value of c. The values of a and b were then chosen to give a reasonable fit to the ex- 

4 

perimentally determined virial coefficients. This was done by calculating the second 
virial coefficient classically from 


B(T) = 27r?7 f [1 - e" V ^/ kT ] r 2 dr (16) 

* / 0 

using various values of a and b to see what combination produced the best agreement 
between the calculated and experimental B(T). Figures 5 and 6 show B(T) calculated 
from equation (16) and using v(r) from equation (15) and also show experimental data 
(refs. 27 and 35). The Slater -Kirkwood potential gives values of V(T) that are up to 
8 percent too low in the range 273 to 1473 K. It will not be used further here. 

Virial coefficients for helium at high temperatures (>500 K) are very insensitive to 
the attractive part of the potential. Even at lower temperatures (down to approximately 
80 K), the shape and depth of the well are inaccurately determined by a virial coefficient 
fit. Thus the attractive well is not accurately determined by the fit of Yntema and 
Schneider. In equation (15) just as in equation (14), the values of potential in the region 
of the well result from extrapolation of the limiting forms for larger and smaller r. 

The well cannot, in fact, be accurately determined by matching virial coefficients 
using the classical formula (eq. (16)). At temperatures where B(T) is sensitive to the 
well shape and depth, a quantum mechanical calculation of B(T) must be made. Figure 6 
shows an example of the inadequacy of the classical formula. It contains B(T) as calcu- 
lated from equation (16) and quantum mechanical calculations taken from reference 36. 
(The 6-12 potential used for the example is from ref. 37. It does not fit the high- 
temperature coefficients very well. For that reason it is used here only to contrast 
classical and quantum results for B(T). ) The quantum mechanical calculation of B(T) 
is much more difficult and lengthy than the classical, and this apparently has prevented a 
quantum mechanical determination of the potential well. 

^It is worth mentioning that London (ref. 33) and Brueckner and Gammel (ref. 34) 
have erroneously reported the Y. S. potential with the constant b given as 4. 82. This 
value yields virial coefficients that are as much as 10 percent too low (between 273 and 
1473 K) whereas b = 4. 72 gives B(T) to within 4 percent (and for most T within 1 per- 
cent) of the measured values. The erroneous value of b yields a potential that is too at- 
tractive. 


22 


14 


12 





□ 


0 


O Yntema-Schneider data 

□ White, Rubin, Camky, and Johnston data 

Yntema-Schneider S 

Weakened-well l Classical calculation (eq. (16)) 

deBoer and MichaelsJ 


200 400 600 800 1000 1200 

Temperature, K 

Figure 5. - Second virial coefficient for helium. 
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re 6. - Second virial coefficient for various potentials. 


In spite of the inadequacy of the classical formula (eq. (16)), another potential has 
been constructed (by the classical eq. (16)) which fits the experimentally measured coef- 
ficients above 100 K as well as the Y. S. potential does). This new potential has a weaker 
well than that of Y. S. and a different core shape. It is introduced here simply to show 
that another potential can give a classical fit to B(T) and yet yield (in the next section) an 
energy spectrum in much better agreement with the experiment. The new potential will 
be designated as the "weakened- well" potential and is shown in figure 4. The classical 
B(T) calculated from it are shown in figure 5. 

Two more potentials are used in the next section to further illustrate the effects of a 
shallower attractive well. These two potentials are identical with the Y. S. potential for 
v(r) > 0 but for v(r) < 0 are equal to oiVy g (r), where a is chosen as 0. 6 for one 
potential and 0. 8 for the other. 

Therefore, a total of four potentials will be used in the next section. Two of these, 
the Y. S. potential and the weakened-well potential, give good fits of the classically cal- 
culated virial coefficients to the experimental ones. The other two potentials with uni- 
formly reduced wells are simply artificial potentials used to show the effects of gradually 
reducing the attractive well. 


CALCULATED SPECTRA 

The integral equations (eqs. (13)) have been solved numerically for each of the four 
potentials: the Y. S. potential, the weakened-well potential of figure 4, and the two po- 
tentials derived from the Y. S. by reducing the well by factors of 0. 6 and 0. 8. The re- 
sulting spectra are presented in this section beginning with the spectrum from the Y. S. 
potential. The qualitative features of that spectrum will be seen to be correct. The en- 
ergies of all excitations will be seen to be high, however, and simple arguments will in- 
dicate that shallower wells should give better results. The spectrum corresponding to 
the weakened-well potential (which was shown to fit virial coefficient measurements) will 
be seen to be much better but still too high in energy. Lastly, the results for the uni- 
formly reduced wells of 0. 8 and 0. 6 of Y. S. values are given. 

The partial wave components of the reaction matrix were calculated from equa- 
tions (D6) and (D7) as the first step in finding the spectrum. The diagonal elements of 
the first three even-numbered waves, calculated from the Y. S. potential, are shown in 
figure 7. These partial waves were summed through 1 = 4 according to equation (D5b)), 
and then equations (13) were solved for the energy spectrum. Figure 8 gives the result. 
Curves obtained using only one or two partial waves in the reaction matrix are also 
shown. A comparison of the calculated spectrum with the experimentally measured one 
shows the energy scale of the present results to be nearly an order of magnitude too high. 
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Figure 8. - Energy spectrum for Yntema- 
Schneider potential. 

The disparity in scale tends to obscure the important similarities: the phonon-like 
low-momentum excitations and the roton minimum. The momenta at which the roton 
minimum and the relative maximum occur are approximately the same as those of the ex- 
perimental spectrum. Previous attempts to include both singular core and attractive well 
have failed to reproduce even these qualitative features. Brueckner and Sawada’s treat- 
ment (refs. 12 and 13) of the hard core gave a qualitatively good spectrum, and as noted 
in a previous section, gave semiquantitative agreement with the experimental spectrum 
for an appropriate choice of a parameter in the theory . But best agreement was achieved 
for a nonphysical value of the parameter that implied that the number of particles in exci- 
ted states exceeded the total number of particles. The attempt by Parry and ter Haar to 
use Brueckner’s method and to include an attractive well led to the loss of even the quali- 
tative features of the helium II spectrum. In the context of these previous results, the 
qualitative features of the present spectrum - phonons at low momentum and roton 
minimum - are gratifying. 

The speed of ordinary sound (first sound) is equal to the initial slope of the energy 

versus momentum curve, because for low-momentum phonons e = pc where c is the 

2 2 

speed of sound. The plot in figure 8 is of E = 2ma e/K as a function of ka where k is 
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the wave number and a is a reference dimension, taken here to be 2. 3 A (2.3x10 m) . 
In terms of the slope dE/d(ka), evaluated at the origin, 

_ de _ fi 2 a dE _ h dE 
dp 2ma 2 K d (ka) 2ma d(ka) 

The value of c from the Y. S. potential is 2300 meters per second, to be compared with 
the actual value of 240 meters per second extrapolated to T = 0. 

It is easy to see what causes the integral equations for E(ka) to give such high en- 
ergies. It is primarily the influence of jf QO in the equations for F(x) and H(x). jt QO 
is a rather large negative number for the Y. S. potential because that potential is nearly 
attractive enough to produce a zero- energy bound state. A potential that is just strong 
enough to have a zero-energy bound state will have a scattering length of -<*>; and by ap- 
pendix F, will also be -°°. For large negative values of Jt QQf the integrals in equa- 

tions (13) may be neglected and the equations become approximately 

F(x) = x 2 + Sirphjq £ 
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H(x) = ftrP^ 


E(x) = YF 2 (x) - H 2 (x) 

Near x = 0, Jt and jf have the following expansions: 

X. X. ox 

2 2 


^xx =jr oo +dX 
2 2 


•*OX = ^oo + bx 


Hence 


E(ka) = V(F + H) (F - H) 


V- 


-1287T 2 P^(2d + b)jf ka 
o' ' oo 


for small ka. The initial slope of the spectrum is thus approximately proportional to 

But study of figure 9 shows the value of d increases rapidly as decreases. 

(The same is true for b. ) Thus the initial slope of the spectrum is very roughly propor- 
tional to -•^ 00 > for large negative jr . 

Figure 9 shows the extreme sensitivity of jt QO to the strength of the attractive well. 
The curves in that figure correspond to the potentials of figure 10. These potentials are 
identical in the core region but have potential wells of three different strengths. The two 
potentials with reduced well strength were obtained from the Y. S. potential by multiply- 
ing all negative values of v(r) by a parameter a, having values of 0. 8 and 0. 6. The re- 
duction of well strength to 60 percent of the Y. S. strength reduces jf Q0 to about 1/40 of 
the value it has for the Y. S. potential as shown in figure 9. This strong sensitivity of 

to well strength is what prompted the construction of the weakened-well potential of 
figure 4. This potential was designed to have a much smaller jt QQ than that of Y. S. and 
yet to fit the measured second virial coefficients just as well. 



Figure 9. - Dependence of diagonal elements of K-matrix on well 
depth (S-wave only). Well-strength parameter a is defined by 
v(r) = av YS (r) for v(r> = 0. 



L _J l 

1 2 3 
r/a (a = 2. 3 A or 2. 3xl0“ 10 m) 


Figure 10. - Potentials with uniformly 
reduced wells. Potentials are iden- 
tical for v(r)>0. For v(r)<0, 
potentials obey v(r) = av Y $(r). 
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The weakened-well potential gives much better results than the Y. S. potential. Fig- 
ure 11 shows the first three even partial waves of the diagonal elements of the reaction 
matrix. Comparison with figure 7 shows that Jt QQ is reduced to one- sixth of the Y. S. 
value. The consequent improvement in the energy scale of the spectrum is evident in fig- 
ure 12. The improvement is significant but short of what is needed to agree with experi- 
ment. Comparing the spectrum with that of the Y. S. potential reveals a shift toward 
lower momentum of the roton minimum and of the relative maximum. The minimum is 
less pronounced. The fraction of particles in the condensed group of states is 91 percent, 
down slightly from the 93 percent result for the Y. S. potential. 

There is some possibility that effects not taken into account in this work might lead 
to an effective weakening of the potential well. For example, the Hamiltonian upon which 
this work is based includes only two-body interaction terms. It is recognized that nonad- 
ditive three- and many-body interactions (refs. 38 and 39) exist in liquids because of the 
composite nature of the atoms. That is, because atoms are not simple "particles" but 
instead have internal structures and because the interparticle force is a result of a mod- 



Figure 11. - Diagonal elements of reaction matrix 
for weakened-well potential. 



Figure 12. - Comparison of spectra from Yntema- 
Schneider and weakened-well potentials with 
experimental spectrum. 
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ification (polarization) of that structure, the force between a pair of atoms is not inde- 
pendent of the presence or absence of other atoms in the vicinity. At low densities this 
is unimportant. But at liquid helium density the effects may not be negligible. Inclusion 
of many-body interactions in the Hamiltonian ^ is out of the question in the present theory, 
but it might be possible to include the many- particle effects approximately by modifying 
the two-particle potential to make it an ’’effective two- particle potential" appropriate for 
the observed liquid density. Whether the presence of the other particles weakens or 
strengthens the Van der Waals attraction between two particles is not at all obvious. 

An approximate microscopic treatment of this problem (ref. 39), which yields a very 
small weakening of the attraction, is discussed in appendix G. Macroscopic methods tak- 



Figure 13. - Comparison of spectra for reduced wells. 
Well-strength parameter a is defined by 
v(r) = av YS (r) for v(r) * 0. 


5 

The distinction between the many-particle interactions under discussion here and 
"many-body terms” in a perturbation expansion must not be forgotten. Even if one could 
exactly diagonalize the Hamiltonian (eq. (2)), one would have in no way included many- 
particle interactions since equation (2) included only two-particle interactions. 
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ing frequency dependence of the dielectric constant into account (ref. 40) have not been ap- 
plied to this specific problem. In any case, because many-body forces or some other 
phenomenon might effectively weaken the attractive well, it is desirable to calculate 
spectra for a series of wells of varying strength. Such an approach serves to uncover 
trends in the spectra and to further show the sensitivity to well strength. 

The potentials of figure 10 form a series of three such potentials, related by the 
single parameter a. The calculated energy spectra can be compared in figure 13. The 
results show that as a decreases the energy scale of the spectrum improves and can 
even fall partially below the observed spectrum. A severe shift toward lower momentum 
occurs, however, which was noted in lesser degree in figure 12. The cause of this shift 
is not known. 

DISCUSSION 

Relation to Previous Work 

To place the results in proper context, the following comparisons of the present 
methods and results are made with other work. 

Most of the microscopic theories have been forced to deal with truncated Hamilton- 
ians. Of these theories several are valid only for weak potentials. Bogoliubov (ref. 10) 
first obtained a phonon- like low-momentum spectrum for weak repulsive interactions 
near absolute zero. His Hamiltonian contained forward, exchange, and pair scattering 
terms, but was diagonalized only by approximating some of the operators by c-numbers. 
Wentzel (ref. 14) and Luban (ref. 15) found an energy gap at low momenta using the pair 
Hamiltonian. The helium II spectrum, of course, has no such gap. Using the idea of 
Girardeau (refs. 22 and 23) that, for an attractive interaction, condensation should oc- 
cur into many states instead of into just one. Luban showed that the spectrum is phonon- 
like at low momentum for an appropriate attractive potential with pair scattering in the 
Hamiltonian. Hence the present work has assumed a generalized condensation and has 
included pair-to-pair scattering. The spectrum obtained herein has a phonon-like low- 
momentum region as observed in liquid helium II. 

The aforementioned works by other authors were based on weak potentials. 

Brueckner and Sawada (refs. 12 and 13) used the reaction matrix method of handling 
strong potentials, but included only forward and exchange scattering and one special type 
of pair scattering in their Hamiltonian. For hard spheres with no attractive well, they 
found a phonon-like low-momentum spectrum. This is similar to Wentzel's result for 
repulsive but weak potentials. The Brueckner and Sawada spectrum had a roton mini- 
mum which approximated the experimental one for an appropriate choice of a parameter 
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in their theory. But this parameter was proportional to the density of condensed (zero- 
momentum) particles, and thus was not really arbitrary. In fact the value of the param- 
eter giving the best spectrum leads to the contradiction that the density of excited parti- 
cles is 2. 7 times the total density. Parry and ter Haar (ref. 17) used approximately the 
same method but consistently handled the density of condensed particles. Their attempt 
to include an attractive well in the potential was unsuccessful. All qualitative similarity 
to the experimental spectrum was lost, including the phonon-like part. 

The present work has used a reaction matrix to handle the strong repulsion. It dif- 
fers from the Brueckner and Sawada reaction matrix, however, in that only kinetic en- 
ergy is included in the propagator G. The additional terms included in the unperturbed 
Hamiltonian by Brueckner and Sawada and by Parry and ter Haar are not necessary with 
net attractive forces. The successful inclusion of the attractive potential in this work is 
due to including pair interaction terms and to assuming the generalized condensation. 

The most significant result of this work is that the two important qualitative features 
of the helium potential - strong repulsion, but net attraction - have been treated with 
methods that were able to produce the two important qualitative features of the excitation 
spectrum - phonons and rotons. 


Discussion of Approximations 

As in all other attempts to derive the energy spectrum, simplications and approxi- 
mations have been made to make the problem tractable. The methods of this work have 
produced an energy spectrum of correct qualitative character from a microscopic theory 
using realistic potentials with both singular core and attractive well. It is appropriate 
to recapitulate the simplications, approximations, and omissions since they may be re- 
sponsible for the lack of quantitative agreement with experiment. Unfortunately the most 
important cause of error has not been identified because of the complexity of the integral 
equations. 

The first simplification was to restrict the second-quantized Hamiltonian to two- 
body interaction terms. Three-, four-, and many-body operators have been omitted, 
but they should be included for an exact treatment of helium. It was shown in a pre- 
vious section that a weakening of the attractive part of the potential by about 40 percent 
gives approximately the right speed of sound. Perhaps the interaction of two helium 
atoms could be modified to this degree by the presence of several near neighbors (speak- 
ing microscopically) or (speaking phenomenologically) by the presence of the dielectric 
medium composed of the other atoms. One treatment of many-body forces, cited in ap- 
pendix G, gives only a 1 percent effect; however, it may not be accurate for liquid he- 
lium. 
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The second simplication was to truncate the Hamiltonian. Only forward scattering, 
exchange scattering, and pair scattering terms were retained. Ways of handling more 
terms than these are not known. Actually the use of K in place of V does, in effect, 
include more terms, but the form of the Hamiltonian is unchanged. 

To allow the inclusion of a singular repulsive core in the potential, the V matrix 
elements were replaced by reaction matrix elements. This, in effect, presummed 
enough terms of the many-body perturbation series to give finite matrix elements in the 
interaction part of the Hamiltonian. But it is shown in appendix B that the perturbation 
series for the free energy derived from the ” Hamiltonian” with K contains some dupli- 
cated terms. The extra terms needed to treat singular cores therefore came at the price 
of including some terms twice. This was recognized by previous authors (refs. 12, 13, 
and 17) but neglected on the basis of canceling errors for ground and excited states. But 
it cannot be said that the elementary excitation spectrum would be unaffected. The dupli- 
cation of terms in the perturbation series for the free energy means the partition function 
would also be in error since 


F = -kT log Z 

The justification for putting K into the Hamiltonian was valid only near absolute 
zero. This is, however, simply a limit on the allowed temperature rather than an ap- 
proximation. If the method is used for elevated temperatures, then an approximation is 
thereby made. 

The Wentzel TEH method is not an approximate method in itself. The equivalent 
Hamiltonian has exactly the same partition function as the Hamiltonian from which it was 
derived and thus the same thermodynamics. It has been argued in an earlier section that 
the spectrum of elementary excitations is also the same if it turns out to be temperature 
independent. 

The last approximation made was to replace the two-body reaction matrix elements 
with their approximately equal one-body central-force counterparts. There is no reason, 
in principle, why this must be done. However, the solution of the integral equation for 
the one-body K matrix was barely practical. The numerical solution of the integral 
equation for a two-body K, which would be a function of four variables instead of two, 
would be impossible without completely different techniques. 

The assumed generalized or "smeared” Bose-Einstein condensation is probably not 
an approximation. It has been shown (refs. 22, 23, and 28) to be a consequence of a pre- 
dominantly attractive interaction, which permits the system to lower its energy by 
spreading out the condensate over many zero-energy (degenerate) but distinguishable 
states. 
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SUMMARY OF RESULTS 


The central result of this work is that the qualitative features of the helium excitation 
spectrum have been derived microscopically from an interparticle potential containing a 
strong repulsive core and an attractive well. The strong repulsive core was handled by 
using reaction matrix elements instead of interparticle potential matrix elements in the 
Hamiltonian. The attractive well, which is strong enough in helium to give the potential 
a negative scattering length, was successfully included by assuming that a smeared Bose- 
Einstein condensation occurs. The pair Hamiltonian (which includes forward, exchange, 
and pair scattering interaction terms) was, in effect, diagonalized by the Thermodynami- 
cally Equivalent Hamiltonian method. The equations derived from these methods were 
solved numerically and yielded spectra with linear behavior at low momenta and roton re- 
gions at high momenta. These correct qualitative features have not previously been de- 
rived from a realistic interparticle potential, containing both a repulsive core and an at- 
tractive well. For the Yntema-Schneider potential for helium the calculated spectrum is 
in poor quantitative agreement with experiment. Weaker potential wells were shown to 
improve the energy scale of the spectrum substantially. A potential was exhibited which 
fits second virial coefficients as well as the Yntema-Schneider potential does, but which 
gave a much improved energy scale. Spectra calculated for a series of three potential 
wells of decreasing strengths showed improving energy scale but a shifting of the roton 
minimum toward lower momentum. 

Lewis Research Center, 

National Aeronautics and Space Administration, 

Cleveland, Ohio, April 30, 1968, 

129-02-05-16-22. 
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APPENDIX A 


SYMBOLS 


a reference dimension, chosen 

herein as approximate core 

o n 

size: 2. 3 A (2. 3x10 m) 

a k’ a k creation and annihilation opera- 
tors for plane-wave states 
with wave vector k 

B k , B k auxiliary operators used in 
finding TEH 

B(T) second virial coefficient 

C+ C k auxiliary operators used in 
finding TEH 

c speed of first (ordinary) sound 

E(k) with argument k: energy of 

state with wave number k 

E(x) with argument x or y: dimen- 

sionless energy E(x) = 
e k x 2ma 2 /K 2 

F Helmholtz free energy 

F Helmholtz free energy for non- 

interacting particles 

f k auxiliary function in integral 

equations for spectrum 

G Green's function 

H general second-quantized Ham- 

iltonian including two-body in- 
teractions 

H q unperturbed Hamiltonian 

Hp pair Hamiltonian 

Hrpg thermodynamically equivalent 

Hamiltonian 


Hj that part of Hp that does not 

contribute to thermodynamics 

H T perturbation Hamiltonian 

h k auxiliary function in integral 

equation for spectrum 

K Planck's constant/277 

K reaction matrix or operator 

Jf dimensionless reaction matrix, 

K X 2mf2/87rfi 2 a 

k wave number; wave vector with 

vector sign deleted for simplic- 
ity 

l angular momentum quantum num- 

ber 

M n n 1 semiinvariant 

m mass of helium atom (or other 

boson) 

N number of particles in system 

N q number of particles in the con- 

densate 

n Q Avogadro's number 

P number of particles in cube of 

side a, P = Na^/fi 

P Q number of condensed particles in 

cube of side a, P Q = N Q a 2 /f2 

p (without subscript) momentum 

p' momentum at which roton mini- 

mum occurs 
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Pi 

wave number; sometimes momen- 

X 

dimensionless wave number, x = ka 


tum 

Z 

canonical partition function 

Po 

nonzero wave number below which 
smeared condensate forms 

Of 

well- strength parameter, a = 1. 0 
for Y. S, potential 

2 

grand partition function 

/3 

1/kT where k is Boltzmann’s 

r 

radial coordinate in spherical co- 


constant 


ordinate system 

y 

r/a 

T 

absolute temperature 

e k 

energy of elementary excitation 

U 

c- number appearing in TEH 


with wave number k 

U o 

c-number appearing in TEH 


variational parameter, finally set 

V 

two-particle interaction potential 


equal to ( a^a£) 


energy operator; one-particle 


chemical potential 


potential energy operator 

^k 

variational parameter, finally set 

V o 

potential inside strong uniform 


equal to ( a j^ a k) 


repulsive core 

O 0 

P 

density, N/fi 

nr) 

v(r) x ma /K 

n 

volume of system 

v(r) 

v(q) 

spherically symmetric interaction 
potential energy function 

three-dimensional Fourier trans- 
form of v(r) divided by £2, for 
momentum transfer q 

w k 

energy of free particle with wave 
number k 
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APPENDIX B 


USE OF K INSTEAD OF V 

Consider the expansion of the free energy F in a series of "semiinvariants. " (See 
for example, refs. 41 and 42. ) Let the exact nontruncated Hamiltonian be split into the 
unperturbed part H Q (kinetic energy) and the perturbation H' (interaction terms), that is, 

H = H q + H’ 


Then F - F Q = M n /n! where F q is the unperturbed free energy and 
n=l 


1 W* 

M< = — / <H'(X))dX 

1 jS Jo 


M 2 = - dx f X dX'[(H’(X)H’(X')> - <H’(X)> (H'(X'))] 
* /3 So So 


Mo = — (* dX f X dX' f X dX"[(H , (X)H , (X')H'(X")) - 3(H , (X)H'(X')> (H’(X")) 
0 0 J) JO JO 


+ 2(H'(X)> <H*(X’))(H'(X")>] 


The bracket means thermodynamic expectation value and is defined as follows: 



where "tr" may be taken to mean a sum over all sets of occupation numbers of the eigen- 
states of H q . H’(X) is a temperature analogue of the interaction picture of an operator. 
Its definition is 
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I 


AH 

H’(A) a e ° 


AH 

H’ e 


The operator H’ for a complete (nontruncated) second-quantized Hamiltonian including 
only two-body interactions is 


H ' = 2 X <k l k 2 |V|k 3 k 4> a k 1 a k 2 a k 3 a k 4 < B1 > 

k l k 2 k 3 k 4 

, T 

as in equation (2). The a s and a's, respectively, create and annihilate plane-wave 
states with the indicated wave vectors. 

The expectation values ( C ) are exceedingly complicated for T i= 0 (finite /3). Also 
the replacement of V by K in the Hamiltonian is valid only at or very near zero temper- 
ature. Therefore, only the limit of zero temperature (j3 = (1/kT) — °°) is considered in 
this work. The excitation spectrum, which is the primary result of this work, is thus 
valid only near T = 0. The use here of the TEH method of Wentzel, is justified, and in 
fact the meaning of ’’excitation spectrum" is definite, only if the spectrum is tempera- 
ture independent for a range of temperature near T = 0. These points are discussed in 
the body of this report. 

At arbitrary temperature 


{€) = 



\Ce 


"$H 


n 0> n l’ n 2’ 


^ ^ (v n r n 2 j ■ ■ ■ 



n 0 » n i> n 2’ 


) 


where ^ indicates a sum over all possible sets of n^, n^, ng, . . ., and 
{n,} 

i = 0, 1, 2, . . . indicate various free-particle states, that is, eigenstates of H Q , the 

kinetic energy. The zero subscript here indicates the zero-momentum state. As T — Q 

-0H O 

and hence 0 — °° the factor e in the numerator causes all other terms to become 

negligible compared to the ground state term, nQ = N and m = 0 for i > 0. Thus 

f* 

Note there is a clear distinction between two-body interactions and two-, three-, 
or n-body terms in a perturbation expansion. 
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lim (G) = lim (G) = — — - ' - I (L N > °> °> • •) = /pi# 

T— 0 j3 — 00 (N, 0, 0, ... 1 1 |N, 0, 0, . . .) 

where 1 0) = | N, 0, 0, 0, . . . ) denotes the ground state of a system of noninteracting 
bosons, which is all particles in the zero-momentum state. 

The first few semiinvariants can now be evaluated in the zero-temperature l im it 




/ XH n j. . -AH 

dA (° | e a^a^e 1 0> (k^ | V |k 3 k 4 ) 


y 

2 8 — 00 8 / y 

k l k 2 k 3 k 4 


/ 


dA 


<°l' 


‘Kj+n 


■w,. 


■O),. 


4 a k 1 a ic 2 a k 3 a k 4 l°^< k l k 2 I 


Vlkgk^ 


U -4 I f 

lr k" k" k *^0 


k 1 k 2 k 3 k 4 


Here o>, is the kinetic energy of a particle with wave number k.. But 

i 1 

^0 I a k ^ \ I = 0 except for kj = kg = k 3 = k 4 = 0. Further o>q = 0 and 

(°l a 0 a 0 a 0 a 0 l°> = <°! a J( a 0 a J " 1 ) a 0 l 0 ^ = N q “ N 0 = n2 " N ~ n2 for N lar S e * Then 


lim M 1 = i lim -i dA N 2 (0| v|0> = - nNv(0) 
3—oo 1 2 8—0° 8 SO 2 


2 j3— oo 8 SO 


(B2) 
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where v(0)/S2 = (0|v|0). The second semiinvariant becomes, as j3 — °° 


lim M, = - lim - dX /** dX»[(H , ^)H , ^ , )> - (H’(X)> (H’ft’))] 

3-»oo /3 -m3 *'0 



dA' 



( k l k 2 I V I k 3 k 4 > < k 5 k 6 I V I k 7 k 8 > 


k l k 2 k 3 k 4 

k 5 k 6 k 7 k 8 


X CO, +CO, -CO. -CO. X f [CO, +C0 -CO. -CO, 

. k ! k 2 k 3 kj \ k 5 k 6 k 7 kg 
x e e 


x 


(° I I °) - (° I l°X° I °> 



x (°l a o a o a k,\. a k/k,, a o a ol 0 ) - {° 0 l v l 00 ) 2 < 0 ! a o a o a o a ol°> 

O 4 D D 


where only kj = k 2 = k^ = k g = 0 appear since other terms have zero expectation value 
in the noninteracting ground state. Only those elements of the V matrix which conserve 
total momentum give nonzero contributions, so 
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i 


J 




X 


< 0 1 a o a O a -k a k a -k a k a O a O I °> 


O' 


' \ u I a o 0 




i 11m i /*\ C dx(i- V N(N . X) 

4 /3— °° 0 / / o. 2 / y 

*^0 •'O L k^O J 


1 p(N ~ 1 ) ii m 1 v 2 (k) 

2 ft 0— 00 0 


k^O 



/ 


-XR 2 k 2 /m/ e Xfi2k2/,m - 1 


dA e 


fi 2 k 2 /m 


■ flN lim i 
2ft 0-*-°° 0 



mv 2 (k) / f , e-^ 2 / 1 ” - 1 


* 2 , 2 
n k 


fi 2 k 2 /m 


k #0 


.flN 

2 


JL A 
8 * 3 / 


3 k mv 2 (k) 
R 2 k 2 


= .£NJ_ m /'° 0 ^2 (k)dk 

2 2 7T 2 K 2 7o 


(B3) 
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In appendix H it is shown that 


lim 

8 < 2„) 6 K 4 





, vOej) v(k 2 ) 


d k, 


v(kj 


k 9 ) 


+ 



8 


— 5(0) 
(277) 3 K 4 



(B4) 


The first three terms in the expansion of the free energy of a nontruncated Hamilto- 
nian are sufficient to indicate what type of "reaction” matrix may be used in the truncated 
Hamiltonian Hp to introduce more terms into the expansion. 

Consider first the series expansion of a general element of the reaction matrix with 
respect to two-particle (noninteracting) plane-wave states 

(kjkg | K | k 3 k 4 ) = (kjkg | (V + VGV + VGVGV + . . . ) | kgk 4 ) 

= < k i k 2 l v l k 3 k 4> + ^ (k 1 k 2 |v|k^ 6 ><k 5 k 6 |G|k 7 k 8 )(k 7 k 8 |v|k 3 k 4 > 

k 5 k 6 k 7 k 8 


+ ^ ^ k l k 2 l^l k 5 k 6^ k 5 k 6 I ^ l k 7 k 8^ k 7 k 8 l^l k 9 k 10^ 

k 5 k 6 k 7 k 8 
k 9 k 10 k ll k 12 

x ^ k 9 k io I G l k ll k 12^^ k ll k 12 f V l k 3 k 4^ + • • ■ ( B5 ) 

Many of the terms vanish because V has nonzero elements only between momentum- 
conserving states. 

The plan is to insert the elements of K into the truncated Hamiltonian (eq. (3)) in 
place of the corresponding elements of V. If the free energy is now expanded in a per- 
turbation series with K’s instead of V's in the interaction term of Hp, the series will 
contain in first order many terms which previously appeared in higher order. Each 
higher order of the series with K in the truncated Hamiltonian has terms which the se- 
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ries with V had only in higher orders or not at all. It is therefore hoped that the series 
based on K (in the truncated Hamiltonian) more nearly resembles the exact perturbation 
series based on a nontruncated Hamiltonian. This can be achieved if G is chosen so its 
matrix elements contain the same energy denominators as those which appear in the exact 
perturbation series. The first semiinvariant containing K in place of V will now be ex- 
amined to show what G must be chosen. It will be seen that "zero- energy” propagators 
should be used. 

Clearly the first term in the expansion (eq. (B5)), which is simply the potential itself, 
will make the same contribution to the free energy as the potential would have made. 
Consider, however, the "second-order" terms of equation (B5), 


Y < k i k 2 I V I k 5 k 6 >< k 5 k 6 I G I k 7 k 8 > < k 7 k 8 1 V I ! k 3 k 4 > ’ 

k 5 k 6 k 7 k 8 


and their contribution to M^. 


Note first that only = kg = kg = = 0 make any contribution at all to Mj. Further, 

since momentum-nonconserving elements of V vanish, only terms with kg = -kg and 
= -kg appear. Further, where G is defined using H Q as the kinetic energy alone, 
only terms with kg = k^ will be nonzero. Hence the second-order terms in equation (B5) 
give in Mj (as /3 — °°) 


An 2 ^T (00 | v|k - k)(k - k | G |k - k) (k - k|v|00> =-p- 

2 k 2 0 



E - K 2 k 2 /m 


= -£NJ_ / k 2 v 2 (k) dk ^ 

2 2tt 2 / -E + R 2 k 2 /m 

* / 0 

This term has the same form as the second-order terms in the perturbation series based 
on V if E is set equal to zero. The resulting zero-energy propagator contains the 
same energy denominators which arise in Mg upon integrating with respect to A' (see 
algebra leading to eq. (B3)). In a similar way the third-order terms in equation (B5) 
produce terms in with the right type of energy denominators if E is always set 
equal to zero. This value for E was also used by Brueckner and Sawada (refs. 12 
and 13). 
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A Green’s function is not uniquely specified until the manner of handling the singular- 
ities is given. If the principal value is taken, then the solution of equation (4) is called 
the reaction matrix. If the integration follows a contour in the complex plane which goes 
above one singularity and below the other (for E # 0), then the solution of equation (4) is 
called the T matrix and may be either "T + ” or ”T_, ” depending on which of the two 
possible ways of going over and under are chosen. The principal value should be used be- 
cause the integrals (such as eq. (B6)) which arise from the expansion of K in terms of 
V actually come from sums over intermediate states. These sums should be like those 
that occur in the semiinvariants Mj. The sums occurring there appeared first as the 
sums in the interaction parts of the Hamiltonian (eq. (Bl)). There E^, Eg, ^3* anc * E^ 
are the allowed wave vectors of plane waves. As the thermodynamic limit is taken, the 
equally spaced allowed values of x-, y-, and z-components of any k^ become more and 
more closely spaced. The sum over the states (see eq. (B5)), second term) therefore 
approaches the principal value of the integral (by the definition of the latter). 

As previously noted (refs. 12, 13, and 17), the use of a reaction matrix is seen to 
cause a duplication of some terms in the perturbation expansion. The "first-order ” term 
of equation (B5), which is simply the potential, produces the result (eq. (B3)) in Mg, that 
is, 


-pN 

2 


--f 

2 Jl SO 


mv (k) dk 


2tt K 


But the same contribution is produced by the second terms of equation (B5) in as 
shown in equation (B6). Similar duplications occur in higher order, though not every 
term is duplicated. According to Parry and ter Harr (ref. 17), it can be shown that, if 
the linked cluster expansion is valid, the same error occurs in the energy of an excited 
state as in the ground state energy. Therefore, the total energy of the system is ex- 
pected to be incorrectly given by this model, but if a small amount of thermal excitation 
is allowed, the difference between the excited and ground state energies will not contain 
the error. 
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APPENDIX C 


CALCULATION OF REACTION MATRIX ELEMENTS 


The types of reaction matrix elements needed are 

(p- p|k| p - 
(p- p|k|- PP ) 
(kp | K | kp> 

(kp | K [ pk) . 


(Cl) 


(As usual the quantities k and p are vectors, but the vector symbol is suppressed for 
simplicity. ) The required elements include those in which k or p may be zero. How 
these elements may be calculated will be discussed here. To reduce the calculation to 
manageable size in terms of machine memory storage and computing time, an approxi- 
mation called the center-of-mass approximation is made. It reduces the elements 
(eq. (Cl)) to elements of a one-particle K operator with respect to one-particle states. 

The K matrix needed is defined by the integral equation 


K = V + VGK = V + VGV + VGVGV + . . . 


where V is the interparticle potential. G has been chosen in this work to be a zero- 
energy Green's function 



where H q is the unperturbed Hamiltonian, taken to be simply the kinetic energy alone. 

The elements of K with respect to two-body free-particle states are required, that 
is, 


<k 1 k 2 |K|k 3 k4) 


ik„-r 1 

where |kgk^) is a state with the wave function (1/S2)e 0 
of the K matrix then satisfies the following equation: 


* 4 r 2 


A general element 
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<k 1 k 2 |K|k 3 k 4 ) = <k 1 k 2 |v|k 3 k 4 >+ £ ( k i k 2 l v l k 5 k 6>< k 5 k 6l G l k 7 k 8>< k 7 k 8l K I k 3 k 4 > 

k 5 k 6 k 7 k 8 

The operator G is diagonal with respect to free- particle eigenfunctions because H Q is 
diagonal. Thus, 


< k 5 k 6 l G l k 7 k 8 > = 


k 

k 5 k 7 

H 


5k 6 k 8 

o 


\ k \ k 
K 5 7 6 8 

fi2 ( k 5 +k e)/ 2m 


(C2) 


Hence 


( k l k 2 I K I k 3 k 4 > = < k l k 2 I V l k 3 k 4> " 


y <v 2 ivi k 5 k 6 > 

k 5 k 6 


f^< k 5 k 6 l K l k 3 k 4> 

k l + k 6 


The V matrix elements may be reduced to one-body matrix elements of a single- particle 
central potential 


(kikgMkgkj) = / d 3 r x 


r 2 r 3 r 4> 


3 3 i ( k r r 1 +k 2‘ r 3 _k 3' r r k 4‘ r 2 ) -s 

d 3 rj d' i r 2 e 1 1 * 6 6 v(r x - r 2 ) 


1_ 

S2 


/ 3 i(k l‘ 
dre 


ko)-r 


V(r)6k l +k 2’ k 3 +k 4 


where r = rj - r 2 and the expression has been integrated with respect to the center-of- 
mass coordinate R. The central potential v(r) has the following one-body matrix ele- 
ments: 
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This is the same as (kjkg | v|kgk 4 ) . The matrix elements of G, however, do not thus 
reduce to one-body elements. Let the definition of the center-of-mass (one-body) G be 

G _ 1 _ 1 
’ -H 0 " -K 2 k 2 /m 


where m occurs rather than 2m because the reduced mass is m/2. Matrix elements 
are taken with respect to single- particle plane-wave states. Then 

_ _ \- k 2’ k 3' k 4 = _ 5k r k 2> k 3~ k 4 

/ K 2 (k x -k 2 ) 2 /4m K 2 (k 2 + kg - 

This quantity approaches the form of the two-body element (eq. (C2)) if kj « -kg. It 
equals the two-body element for a pair excitation, where kj and kg are equal and op- 
posite momenta, produced, for example, by pair-to-pair scattering from the condensate. 
For kj and kg such that the total momentum is not negligible compared to the relative 
momentum, the approximation is poorer. The exact effect of using this approximation is 
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n 

certainly not known for the present application. Parry and ter Haar (ref. 17) made 
some approximate calculations to estimate errors and drew the qualitative conclusion that 
the error would be small. 

The matrix elements (eqs. (Cl)) thus reduce (to the accuracy of the center-of-mass 
approximation) as follows: 


(p - p!k|p - p> - <p|k|p) a K pp 


(p - p|K|-pp> - <p|K|-p> 3 K p _ p 

(kp)K|kp> - ^-=-2 |k| s I^_p k-p 

2 2 

(kp|K|pk) - |k| 

2 2 


There are no special problems in calculating the first two of the elements. The elements 
K k_p k-p are reall y of the simple form where q = (k - p)/2 (fig. 14). In the in- 

2 2 

tegral equations derived in the text, these elements appear in an integral of the form 


Q 

(2tt) 3 



d 3 P k Uv) 

k-p k-p 
2 2 



Figure 14. - Change of variables. 


^It has been used in related treatments of the imperfect boson gas (refs. 12, 13, 
and 17) . 
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Note that q = - Vk 2 + p 2 - 2kp/i from the law of cosines, where p is the cosine of the 
angle between k and p. Let the polar axis of the p- space be parallel to k, and denote 

K kzP ± klP =K q±q bY ^ Th6n 

2 2 


n 

8tt 3 


/* d 3 p K ?(P)=— /*°° P 2 dp £(p) Z 1 

/ Kl£±lElP 4 tt 2 *ol A 

J 2 2 


P 2 dp $(p) dp. K.(I Vk 2 + p 2 - 2kp/i 


o A* 2 /*q=(k+p)/2 

^ / P 2 dpi;(p) / dq-2-K ± (q) 

7 r 2 *q) ^=[k-i 


-p|/2 


kp 


n r° 

2 , ^ 

77 k 


P dp 


/*(k+p)/2 

£(P) / , . qK (q) dq 

yk-p/2 * 


77 2 k ^ 


P dp |(p) 


^ 


where 


APPENDIX D 


PARTIAL WAVE DECOMPOSITION OF REACTION MATRIX 

As previously stated, the two-body reaction matrix elements will be approximated by 
the one-body, center-of-mass elements. The integral equation for that operator is 

K = V + VGK 

where G = -1/H Q . The desired matrix elements then satisfy the following equation: 
<k|K|p> = <k|v|p> + (S') VGK | p) 

= ( k| v |p )+ _i_ / d 3 q(k|v|q> ^ <q|K|p> (Dl) 

(2tt) 3 J h 2 q 2 

The calculation is facilitated by decomposing the terms into angular momentum par- 
tial waves. Thus a matrix element of an operator C may be written as 

(^Mp) s ^ (2Z + l)(kMp> z P z (£ * P) 

l 

where £ and p are unit vectors parallel to k and p, respectively, P ^(x) is a Legendre 
polynomial, and 


<k|tf|p)jsi (k|/?|p> P z (k • p) d(k • p) 

The Legendre polynomials have the property that 

r l 25. 

/ P,(p)P m (M) dp =— ^2- 

J- 1 1 m 21 + 1 

Choosing a spherical coordinate system with IF parallel to the polar axis and with p 
having an azimuthal angle of zero (see fig. 15), multiplying the terms in equation (Dl) 
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Figure 15. - Visual aid for partial wave 
expansion. 


A A* 


by - P (k • p), and integrating over all values of k • p between -1 and 1 result in 


i J (k|K|p> P z (k • p) d(P • p) = i ^ (k| v|p) P^(k • p) d(k • p) 


— f d(k-p) f < 

/ J 


d 3 q(k | V | q) (q|K|p)P (k • p) (D2) 


h 2 q 2 


Let (k |k [p) = (2 1 + l)P z (k • p)<k|K|p> z and <£|v[p> = £ (21 + l)P^(k • p)<k|v|p) z ; 

l 

then (k|K|p) z = i (k|K|p)P z (k • p) d(k • p) and 

(k | V | p) ^ (k|v|p)P z (k • p) d(k • p). 

Equation (D2) reads, with several substitutions, as follows: 


50 


I 



(k | K | p) £ - <k|v|p) z 


1 

16tt 3 






X ^ {2V + l)P r (k • q)(k| v|q) z , “0_ 
V K 2 q 2 


x ^ (2Z m + l)P z „(q * p)(q|K|p) ZM P z (k • p) (D3) 
V ' 

where <p is the azimuthal angle of q. Figure 15 shows angular variables which are re- 
lated to k, p, and q by cos 0 Q = k • p, cos a> = k • q, and cos v = p • q. The variable 
v can be eliminated by using the addition formula for spherical harmonics (see, for ex- 
ample, Morse and Feshbach (ref. 43)) 


Z” 

V t ,,(cos v) = e m ~ m || P™, (cos e Q )P™ (cos o>) cos (m <p) e Q = 1, e m = 2, m > 1 

m=0 


or 


P Z”(P 


q) 


L 

-z 

m=0 


m 


(Z" - m)! 
(Z ' ’ + m) ! 


p z"^ 


p)Pj?,(k 


q) cos (m cp) 


(D4) 


But in integrating equation (D3) over cp, only the m = 0 term of equation (D4) will sur- 
vive, so equation (D3) becomes, substituting for P^(p • q), 
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(k | K | p) £ = (k | V | p) j 


— / d(k • p) / dq / d(k • q) \ (2V 

1 l y 


+ 1 )( 2 1 "+ 1 ) 


x (k|v|q> zl ^ (q|K|p) zfl P zf (k • q)P z „(k • q)P z „(k ■ q)P z (k • p) 
fi 2 




(k | V | p) ^ — ~r / d(k ■ p) 

87 T 


2{2l' + l)(k|v|q) 7 , 


mf2 


R 


x (q|K|p) it P r (k • p)P 7 (k • p) 


or 


<k|K|p) z -(k|v|p ) 7 




dq<k|v|q) z ^<q|K|p>j 

R 2 


(D5a) 


TRis last equation may be solved for any l, and the entire reaction matrix may then be 
found by 


<k|K|p> =Yj {21 + ' P>< k l K |P> Z 

l 


(D5b) 


The calculations by machine are simpler when the equations are written in dimen- 
sionless form, and the resulting forms are those required in equation (13). Recall that 

2 o 

— k and ma V. Multiplying equation (D5a) by gives 


47 T J2 
R a 


4 ’ K 2 a 


<*l -s=- K l»> -i ttf v|p >, ' -i f v|,) k|p> , 

47rR 2 a l 4 h 2 a 1 2ir 2 J Q R 2 a l 47rR 2 a l 
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Note that 


( k I^T K |p) = “ f d(k • p) P,(k • P) — f d 3 r(k|r)V(r)(r|p) 

K a l 2 J_ x h 2 a J 

OO 1 

= — / d(k • p) P ? (k • p) — f r 2 dr V(r) f d(k • r) V {2V + 1) 

J_i R 2 a '0 J-l vi m 


/ 2tt 

d<? 


and 


(kjr> = i*(2Z + l)P z (k • r) 

l 


Then 


<k|— V |P> = 27r f d ( k • P) p Z ( k ■ P) — A 2 dr E ( 2z ’ + 1 ) p z .(P • k)j r (kr)j z ,(pr)v(r) 
h 2 a l J-i R 2 a J V 

= 47r — /* r 2 dr v(r)j ; (kr)j ; (pr) 
h 2 a ^ * * 

/ °° 2 

r j ; (kay)j ; (pay)y^y) dy y = - 
Li a 

where j^(x) is the spherical Bessel function of the first kind of order l. So 
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(D6) 


<k|jr|p) z = <k|^|p) z - ^ d(qa) <k | q> ^ (q|/|p)j 

where 

(k[-T|q ) t = J Q rr(y)j z (kay)j z (qay) dy (D7) 

Equation (D5a) or the dimensionless versions (eqs. (D6) and (D7)) can be solved 
rather simply with a computer for any potential for which (k| V|p) exists. The simplest 
method of solution is "guess and iterate. " The partial wave components (k|^|p)^ are 
first calculated in the dimensionless form (eq. (D7)), and then they are substituted into 
the right-hand side of equation (D6) along with the first guess for (k|jf|p)^. The result- 
ing left-hand side could serve as the next approximation to (k|jir|p)^, but averaging that 
value with the original guess before starting a new iteration helps reduce or prevent 
"oscillations" about the final answer. The situation is rather analogous to underdamping, 
overdamping, or critical damping in a mechanical oscillator. The relative weighting of 
the old and new values of the elements (k |jf| p)^ determine the damping characteristics. 
Where the elements (k|^|p) 7 are large and the initial guesses for (k|jf|p) were not 

very close, a large relative weight had to be given to the n approximation in compari- 
th 

son to the (n + l) in in order to prevent overshooting the correct answer and oscillating 
about it. 

For a hard core or any singular repulsive core that increases faster than 1/r, 
(k[')^|p)Q (S-wave) does not exist, and higher partial waves may also be divergent. How- 
ever, (k|jf|p)^ does exist for such potentials for any l if it is considered to be the limit 
of the class of solutions of equation (D5a) as a finite repulsive core is made progressively 
stronger. In fact, this is the way the integral equation was solved by machine. A mod- 
estly strong repulsive core was introduced for the first iteration, and the strength was 
increased each iteration until further increases had negligible effect on the resulting re- 
action matrix elements. Figure 16 shows the Y. S. potential with the core cut off at 
V = 300 and 400. These energies correspond to V/k = 684 and 912 K. Values of the 
diagonal elements, (p|jf|p) 0 (the S-wave component), of the reaction matrix differed by 
less than 1 percent for these two cutoffs. Further increase of the cutoff would have a 
still smaller effect on the elements. The higher partial waves are less sensitive to the 
core details because of the angular momentum barrier. 

The matrix elements of V and K are square arrays of numbers. Storing these ar- 
rays for small enough increments between values and large enough ranges for good ac- 
curacy and doing the required numerical operations in a reasonable time on the computer 
was a real problem. For a given potential (such as Y. S. ), about 5 minutes of computing 
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Figure 16. - Core cutoff of helium 
potential. 


time was needed to calculate the K matrix. Probably a more efficient method exists for 
solving the linear integral equation for K ? although the iterative method is undoubtedly 
the most straightforward. 

From (k|K|p> = £ (21 + l)Pj(k • p)(k|K|p) z it is evident that, if (k|K|p) + 

(k|K[-p) occurs, the odd partial wave terms will cancel out because (k|K|p)^ depends 
only on the magnitude of k and p whereas P^(k • p) = ±P^(-k • p), depending on whether 
l is even (upper sign) or odd (lower sign). If either k or p is zero, then (k|K|p)^ = 0 
for l ± 0, so only the S-wave need be calculated. 
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APPENDIX E 


BOGOLIUBOV TRANSFORMATION 

The thermodynamically equivalent Hamiltonian is easily diagonalized by a Bogoliubov 
transformation. The form to be diagonalized is 


H TE " u + ^ [ f k a k a k + \ h k(‘ 


k 


+ + ' 
a, a , + a , a, 
k-k -kky 


This can be put into the form H =U’ where U’ contains no operators, by 


making the following transformations: 


a k = Vk + V-k 


and 


a k " u k®k + v k a -k 


where u^ = u_ k = u* and v fc = v_ k = v£. The transformation is canonical; that is, the 
new operators obey the same commutation rules as the old, if \ ~ = 1 for all k. 

The transformation brings Hrpg to the form 

H TE = u + X fl t("k“k“k + v k“-k"-k) + h k u k v k(“k“k + “-k“-k) 


^kVk + \ \( U k + V k) ( a k“-k + a -k a k) 


This form explicitly shows the division of the Hamiltonian into diagonal and off-diagonal 
parts. The nondiagonal parts may be eliminated by the following choice for the remaining 
condition on u k and v^: 


*kVk + \ h k( u k 
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for all k. By solving the two conditions on and v^, they can be eliminated from 
H te in favor of f k , h k , and e k> Then 




+ 


1 


k 2 e. 


1 


and 


u k v k = - 7 - 

* e k 


with 



The final result is 


h te = u + ; Z (£ k - V + Z £ k“k“k 

* k k 


If e k = 0, the transformation is invalid. 

Expressions for £ k and 77 k can now be found. Luban shows (ref. 15) that the ther- 
modynamic equivalence of H,p E with holds if £ = (a£a k ) and 77 = (a^a^, where 
averages are with respect to a grand ensemble. The set of states created by the o? k are 
used to compute (a^a k )^ These operators create noninteracting quasiparticles, and they 
obey Bose-Einstein commutation relations. Then 


< a k a k> = u k< a k a k> + v k<"-k“-k> + n k v k< a k a -k + a -k a k> 


= u k( a k a k> +v k < 1 + a -k a -k> 

= ( u k + v k)<°'k Q 'k> +v k 

2 2 

Using u k and v£ from previous equations and noting the usual Bose occupation formula 
(for energy e k measured with respect to the chemical potential) result in 
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and (a ^a^.) is evaluated in similar fashion 

a -k a k = ("k fl -k + v k a k> (u k a k + v k«-k> 


= "k^-k^k + v k“k“-k + Vk(°-k ff -k + “k^k) 


< a -k a k> = Yk^-k a -k + a k a k> 


- u k v k (2a^ k + 1) 



58 



APPENDIX F 


RELATION BETWEEN K Q0 AND SCATTERING LENGTH 

It will be shown here that a x q is equal to the scattering length. The one-body 
reaction matrix used in this work has been defined using a zero-energy propagator. A 

O 

diagonal element of the reaction matrix thus satisfies 


<p|k|p> = ( p | v | p) 



( P I V | k) 


m 


-fi 2 k 2 


(k | K | p) 


The more commonly used reaction matrix K* has a nonzero propagator, that is, 


(p|K'|p) = < P | V | p) + ) ( P | V | k> m <k|K»|p) 

L-j k 2 ( p 2 - k 2 ) 

k 

This reaction matrix K' has the property that 

<p|K’|p> = (p|v|^ p ) 

where V is the central potential, and |i^ p ) is an exact state of momentum p in the 
presence of the potential, and which is asymptotic to a plane wave plus a standing radial 
wave. (See, for example, Thaler (ref. 44).) Thus, 

(p|K’|p> = ^ J d 3 r ^ p (r)V(r)^ p (r) 

where <P p (r) is a plane-wave function, ^ p (r) is the exact wave function, and £2 is the 
volume of the box containing the particle. For the limit as p — 0 only the 2 = 0 partial 
wave contributes. Thus 


8 


Vector signs on momenta are suppressed. 
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lim (p|K f |p) = lim — f r 2 dr ^>^(r)V(r)R^(r) 
p— 0 p-*0 f2 J ^ " 

= — lim f— tan 6^) 

to m p— oyp J 

= i2L K? L 
fi m 

where R^(r) is the S-wave radial wave function, cp^(r) is the S-wave radial wave func- 
tion of a plane wave, L s lim | - tan is the scattering length, and 6$ is the S-wave 

P-°° Lp J 

phase shift for momentum p. Now, the diagonal element of K’, (p|K f |p), is continuous 
at p = 0; so 

<0 |K T 1 0) =— — L 
to m 

From the defining equations for K and K’, the elements (OjKjO) and(0|K' |0) can be 
expanded as follows: 


and 


< 0 1 K 1 0> = ( 0 1 V 1 0) 



-fiV 


(k | V 1 0) 


+ . . . 


(0 | K T 1 0) = (0| v|0) 


I 


(0|v|k) 


m 


-R 2 k 2 


(k | V 1 0) 


+ 


Thus (0|K|0) = <0 |K T 1 0). Then 
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<o|jf|o) =— — <o|k|o> 

4 ’k\ 

= lR^ (0lK , j0} 

a 4w ^2 

= L 
a 

Hence, 

a x jt QQ = a( 0 1 Jf 1 0) = L 



APPENDIX G 


APPROXIMATE EVALUATION OF SINANAGLU EFFECTIVE 
TWO-PARTICLE POTENTIAL 

Only two-particle interactions have been included in this work. Sinanaglu (ref. 39) 
has given an approximate way to modify the two-particle potential V(r) to include the in- 
fluence of the other particles in the liquid. He finds that the attractive well of the inter- 
particle potential may be weakened as much as 10 to 40 percent by the medium. The 
largest effects occur for different solvent and solute species. But for two atoms of the 
same element in the liquid of that element, he derives an effective potential of interaction 
between two atoms as 


V eff (r) - 


vacuum 


(r) x B(r) 


where to sufficient accuracy for helium: 


B(r) - 1 - 





where n^ is the zero-frequency index of refraction, a is the position of the zero of 
V(r) if the vacuum potential is represented by the Lennard-Jones 6-12 potential. The 
function L(r,o) is given graphically (ref. 39) from numerical calculations. It is es- 
sentially equal to 2 for (r/cr) > 3. For 1. 6 < (r/cj) < 3, L(r,cr) & — (— - 1.3); for 

9 3 \cr / 

0. 5 < (r/cr) <1.6, L(r, a) & 0. 3[(r/a) - 0. 5] ; and for (r/cr) < 0. 5, L(r, cr) = 0. For 
values taken directly from the graph of reference 39 for L, the function B(r) differs 
from one by about 1 percent at the most. The main reason there is such a small reduc- 
tion of strength for helium (compared with an order of magnitude more in some other 
liquids considered by Sinanaglu) is that n^, the index of refraction, is only 1. 03 for he- 
lium whereas for many liquids it is 1. 3 or 1. 4. Thus n^ - 1 is only one-tenth as large 
for helium as for many liquids. 
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APPENDIX H 


THE THIRD SEMIINVARIANT 

The third semiinvariant has been defined as 

1 /*£ /*^ ' 

M, = — / dX / dX’ / dX ,, [(H'(X)H’(X , )H’(X ,T )) - 3{H , (X)H , (X')> (H'(X ,r )) 

0 j 8 ^0 */0 

+ 2(H'(X)) <H’(X’)> <H’(X")>] 

where 




k l k 2 k 3 k 4 


<kl k 2 l v l k 3 k 4 >< 1 < 2 a k3 a k4 


and 


(kjkglV k 3 k 4 >e 


^%«W“ k 4 , 


ki k 2 k 3 k 4 


Only contributions from the (H'(X)H , (X , )H , (X")} term will survive cancellation; so con- 
sider only that term, which will be denoted by M^, 
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dA"(H'(A)H'(A')H'(A")) 




x(kj - kj | V|k 2 - k 2 ) <k 2 - k 2 Jv| 00 > 


fli 2 /m) 


2Xk j +2X '(k^-k^+2 A ' .k| 


x (°| a O a O a ' 


kj -k 1 a k 1 a -k 1 a k2 a -k 2 a k2 a " k 2 a ° a 


,1°) 


iVi A A- t 

2 k*0 ^ ^ 


dA"(00 1 V|k - k) (kO| v|kO) (k - k|v|00) 


x e^ 1 +2A M k )/n| a + a + a a a + a + a, a a + a + a a l0\ 

e < u I a O a O k a -k a k a (A 0 k a -k a O a 0 1 

+ terms that cancel with other parts of the semiinvariant 


From this point pursue only the terms that do not cancel. These terms become 



v(k,) v(k. - k„) v(k.) (n /m) 

dA” i 1 “ L e 


-2xkj+2A’(kj-k2)+2X"k; 


/-/ 


dv ,y(k) v(0) v(k) e (R 2 /m)(-2Xk 2 +2X"k 2 ) N 3 

n n n 


where terms proportional to powers of N less than one have been dropped and it has 
been recognized that < 0 I a, 1 0 ) =0 for k * 0 . Further manipulation yields 
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v(kj) v(k 2 ) 



v(k. -l) + ‘l5!!!!!5(0) 
1 z 8 n n 4 



k*0 
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